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Abstract 

We propose that a regulation mechanism based on Hebbian covariance plasticity 
may cause the brain to operate near criticality. We analyze the effect of such a 
regulation on the dynamics of a network with excitatory and inhibitory neurons 
and uniform connectivity within and across the two populations. We show that, 
under broad conditions, the system converges to a critical state lying at the common 
boundary of three regions in parameter space; these correspond to three modes of 
behavior: high activity, low activity, oscillation. 
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1 Introduction 



That evolved brains are highly sensitive organs is an everyday observation. Viewed as a 
dynamical system, a brain may be said to be unusually susceptible to perturbations and 
initial conditions. This leads one to ask whether brains may be operating near some form 
of instability, or criticality, a hypothesis related to the notions of computation at the edge 
of chaos (Langton 1990) and self-organized criticality (Bak et al. 1987). In this paper 
we propose that while most regulation mechanisms at work in the brain act according 
to a classical homeostasis schema, i.e., have a stabilizing effect, an opposite effect could 
result from the regulation of synaptic weights by a specific form of Hebbian covariance 
plasticity. Such a regulation may bring the system near criticality. We suggest that 
regulated criticality may be the mechanism whereby sensitivity is maintained throughout 
life in the face of ongoing changes in brain connectivity. 

Hebbian synaptic plasticity (Hebb 1949) plays an important role in the development 
of the nervous system, and is also believed to underlie many instances of learning in 
the adult. A covariance rule of Hebbian plasticity roughly states that the change in the 
efficacy of a given synapse varies in proportion to the covariance between the presynaptic 
and postsynaptic activities. As noted by many authors (e.g. Sejnowski 1977a, 1977b; 
Bienenstock et al. 1982; Linsker 1986; Sejnowski et al. 1988), a covariance-type rule 
is preferable to a rule that uses the mere product of pre- and post-synaptic activities 
because the covariance rule predicts not only weight increases but also activity-related 
weight decreases, and as a consequence allows convergence to non-trivial connectivity 
states. Some forms of covariance plasticity have been shown to be optimal for information 
storage (Willshaw and Dayan 1990; Dayan and Willshaw 1991; Dayan and Sejnowski 
1993). Also, evidence for Hebbian plasticity of the covariance type has been reported in 
several preparations (Fregnac et al. 1988, 1992; Stanton and Sejnowski 1989; Artola et 
al. 1990; Dudek and Bear 1992). 

We shall investigate, in a simple network including excitatory and inhibitory neurons, 
the effect of covariance plasticity acting as a mechanism of regulation, rather than super- 
vised learning. Synaptic modification results in changes — quantitative or qualitative — in 
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the activity that reverberates in the network; these changes in turn cause further modifi- 
cation of the weights, thereby creating a feedback loop between activity and connectivity. 
Studying this loop as such, i.e., independently from any input and output, we demon- 
strate that, under fairly general conditions, it causes the network to converge to a critical 
surface in parameter space, the locus of an abrupt transition between different activity 
modes. In Metzger and Lehmann (1990, 1994) a similar Hebbian rule has been studied 
in the context of supervised learning of temporal sequences. 

Schematically, the convergence to a critical state can be explained as follows. One 
mode of behavior of a network including excitatory and inhibitory neurons is oscilla- 
tion; such behavior takes place if the synaptic weights linking excitatory neurons to each 
other — we will refer to these as E-to-E weights — are high enough but not too high. Os- 
cillation entails high covariance values, hence, according to the covariance rule, results in 
further increase of the E-to-E weights. If however these weights are higher than a certain 
critical value — which depends on other parameters of the system — oscillatory behavior 
is impossible, hence covariance is low or zero, hence, in accord with the covariance rule 
used, the E-to-E weights decrease. As a result, the E-to-E weights stabilize around the 
critical surface that separates the region of oscillation from the region(s) of steady firing. 

Our study is conducted in the simplest type of network that will support oscillatory 
activity: all synaptic weights of a given type are given identical values, and so are all 
firing thresholds of a given type. This results in a system with just six parameters — four 
synaptic weights and two thresholds — and a limited range of behaviors. Essentially, all 
neurons fire uniformly, either at a constant rate (the number of possible rates of firing 
is one or two, depending on parameters) or periodically in time. In the thermodynamic, 
i.e., large-size, limit, the dynamics of the network is adequately described by a system of 
differential equations obtained through a classical mean-field approximation. 

We first perform a simple bifurcation analysis of this differential system (Guckenheimer 
and Holmes 1983). We then show that the effect of covariance regulation is to stabilize 
the parameter state at a surface of transition, where the dynamics exhibits an instability. 
Such a critical parameter state for a dynamical system may be characterized as degenerate. 
A generic, i.e., non-exceptional, state is one where one expects to find the system in the 
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absence of further assumptions. Mathematically, a generic state of a dynamical system is 
in the interior of a parameter region corresponding to a given behavior, and the system in 
such a parameter state is said to be structurally stable; the set of non-generic parameter 
states has measure zero. We shall show that a state of higher degeneracy, characterized 
as a point of intersection of several critical surfaces, can be achieved by the simultaneous 
regulation of several parameters. In the vicinity of that highly degenerate state, the 
system displays a range of behaviors, including chaos. 

The plan of the paper is as follows. In the next section we study the dynamical 
properties of our simple network — in the differential-equation formulation — with fixed 
parameters (synaptic weights and firing thresholds). We characterize the bifurcations 
which take place at the boundaries between domains corresponding to different modes 
of behavior. This study is conducted for a reduced system, where the thresholds are 
eliminated in such a way as to render the dynamics symmetric about the origin. Section 
^| describes the regulation equations. Section |] describes the behavior of the regulated 
reduced system, and Section |5| that of the regulated full system. 

2 The fixed-parameter model 

This section describes the dynamics of the model with fixed parameters. We first briefly 
describe a network consisting of a large number (2N) of binary- valued neurons operating 
under a stochastic dynamics. However, rather than using this network for our study of 
plasticity, we make a number of simplifications and approximations, leading to a deter- 
ministic two- variable differential system with just six parameters. The two variables are 
the excitatory and inhibitory population averages of cell activity in the 2iV-dimensional 
model; the six parameters include the four average weights of the synapses within and 
between these two populations, as well as the average firing thresholds for the two pop- 
ulations. We then study the asymptotic behavior of this differential system for various 
parameter values. Different types of asymptotic behavior, in different regions of the 
parameter space, correspond to different phases of the stochastic system, and we pay par- 
ticular attention to the bifurcations of the solutions, where the bifurcation parameters are 
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the synaptic weights — see Schuster and Wagner (1990) and Borisyuk and Kirillov (1992) 
for a related bifurcation analysis. Bifurcations correspond to phase transitions in the 
statistical-physics formulation (the original 2iV-dimensional model). 

We consider a fully-connected network of TV excitatory and N inhibitory linear-sigmoidal 
{0, l}-valued neurons,[] operating under a stochastic dynamics. We denote the activity of 
the i-th excitatory, resp. inhibitory, neuron by xf (t), resp. x\(t), with xf(t),x\(t) G {0, 1}, 
i = 1 . . . N, and we denote the synaptic weights by wf E , wfl, wjf, wjj, i,j — l...N, where 
i is postsynaptic and j presynaptic, and the superscripts indicate, for each of the two neu- 
rons, whether it is excitatory or inhibitory. Thus, for all % and j, wf E and w lE are positive 
or zero, whereas w E } and w)] are negative or zero. 

The local field on excitatory neuron i, i.e., the difference between its membrane po- 
tential and its firing threshold hf, is gf (t) = J2j w T^ x f(t) + Hj w fj x )if) — hf. Similarly, 
the local field on inhibitory neuron j is g\{t) = J2j w }f x J(t) + J2j w jj x }('t) ~ h\, where h\ 
is the threshold of inhibitory neuron i. The network dynamics is defined by: (i) select- 
ing at random, with uniform probability, one of the 27V" neurons; (ii) computing its local 
field g(t), of the form gf(t) or g\{t); and (iii) defining the state of the network at time 
t + 5t to be equal to the state at time t except, possibly, for the selected neuron, whose 
state becomes — or stays — 1 with probability |(1 + tanh(/3g(£))). Parameter (3 is a fixed 
non-negative number, an inverse temperature. The temperature T = 1/(3 measures the 
amount of noise in the system: the higher the temperature, the noisier the dynamics. The 
update interval is 5t — 1/ (2 AT), so that each neuron is updated on average once every time 
unit. This asynchronous dynamics, of the Glauber type (Glauber 1963), is widely used 
in statistical-mechanics models; it lends itself to a convenient mean-field approximation 
(see below). 

A system such as the one just described will exhibit a highly diverse range of behaviors, 
depending on the values of the synaptic weights and firing thresholds. But we now make 
the much simplifying assumption that synaptic weights and firing thresholds are uniform 
across each class. Specifically, for all i, j = 1, . . . , N, we assume that hf = h E , h] = h l , 
u>ee = W EE /N, wf} = -w EI /N, wjf = w 1E /N, and wj] = -w u /N, where h E , h 1 , w EE , 
1 It is not essential that the numbers of excitatory and inhibitory neurons be the same. 
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w , w and w are fixed parameters, and w , w , w and w are non-negative. The 
dynamics is thus parameterized by six constants, four synaptic weights and two thresholds; 
(3 is a mere multiplicative factor common to all six parameters, yet it is convenient to use 
it as a seventh parameter. Unless otherwise mentioned, (3 will be 1. 

Due to this uniformity assumption, all neurons in any of the two populations experience 
the same field at any given time. This system exhibits a limited number of fairly simple 
behaviors, of which Figure 1 is an example. This figure shows the time variation of (xf(t)) 
and (x\(t)), the average activation levels across the excitatory and inhibitory populations. 
In this example, parameters are: N = 70, w EE = 12, w 1E = 8, w El = 10, w 11 = 2, 
h E = 1, h — 3. One unit on the time axis corresponds to 2N updates, so that each 
neuron is updated, on average, once every time unit. For these parameter values, the 
system oscillates. Note that the oscillation is not perfectly regular, a finite-size effect. 
Note also that the inhibitory activity lags somewhat behind the excitatory activity: the 
excitatory neurons first trigger the inhibitory ones, which in turn extinguish, for a while, 
the excitatory population. 

(Insert Figure 1 around here) 

The presence of oscillations and the amplitude and shape of the waveform depend 
on the various parameters. However, rather than pursuing this study of the stochastic 
system, we shall consider the approximation that obtains in the thermodynamic limit, 
that is, when iV — > oo. The update interval 5t = 1/(2N) then goes to and so does 
each individual synaptic weight. Straightforward approximations (Rubin 1988; Schuster 
and Wagner 1990) then lead to a continuous-time differential system for the population 
averages of the excitatory and inhibitory activation levels, which we denote, respectively, 
by s and a: 

s(t) = .5 - s(t) + .5tanh[/3(u> EE s(t) - w EI a(t) - h E )} 

(!) 

a(t) = .5 - a{t) + .5 tanh[/?(w IE s(t) - w ll a{t) - h 1 )}. 

Note that the variables s(t) and a(t) remain at all t within the interval [0,1]. When (3 = 
system |I| has a unique attractor, (s,a) = (.5, .5). Indeed, in the high-temperature limit, 
all neurons act independently of each other and fire with probability .5 at each time. 
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We shall now make a last simplification, whose purpose it is to render (.5, .5) a fixed 
point — though not necessarily stable — at all temperatures and for all values of the synaptic 
weights. This is easily achieved by letting the thresholds h E and h l be determined by the 
synaptic weights as follows: 

h E = .5(w EE - w El ) 

(2) 

h l = .5(w 1E - w 11 ). 

It is then convenient to adopt the change of variables: s i— > s — .5, a i— > a — .5, and system 
D becomes: 

J s(t) = -s(t) + .5tanh[/3(w EE s(t) - w EI a(t))} 
I a(t) = -a(t) + .5 tanh[(3(w 1E s(t) - w ll a(t))} . 
In 0, the variables s and a are in the interval [—.5, +.5], and the only parameters left are 
the four synaptic weights and the inverse temperature. For all parameter values, the origin 
is a fixed point of system |3|. A different position for the fixed point could be obtained with 
an appropriate modification of equations |2], yet in the current version the fixed point is 
also a center of symmetry. For the moment, this hard-wired symmetry should be regarded 
as an ad-hoc device, whose purpose is to make the analysis more convenient. We shall 
refer to system [1] as the full system, and to system ^ as the reduced system. We shall see 
in Section [| that under appropriate regulation the two systems behave similarly. 

We now discuss some important properties of the reduced system, system |3| (see also 
Rubin 1988). Consider first Figure 2a (phase diagram), which shows four trajectories of 
the state (s(t),a(t)); the starting points of these trajectories are indicated by triangles. 
The parameters (synaptic weights) used in this example are identical to those used in 
Figure 1, i.e., w EE = 12, w lE = 8, w El = 10, w 11 = 2. As expected, the asymptotic 
behavior is periodic; there is a limit cycle which attracts all points of the square [—.5, .5] 2 , 
except the unstable equilibrium (0,0). Motion is counterclockwise, for, as mentioned 
above, a(t) lags behind s(t). 

In addition to these four orbits, Figure 2a shows two curves, the s- and a-nullclines 
for system [| These are the loci of the points (s,a) such that ds/dt, resp. da/dt, vanish. 
The equations for the s- and cr-nullclines are easily seen to be, respectively: 

a = ^-(w EE s-T tanh' 1 (2s)), (4) 
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s = ^-(w u a + Ttanh~ 1 (2a)). (5) 
w 1E 

The cr-nullcline is an increasing sigmoid-shaped curve, whereas the s-nullcline generally 
has the shape of an 'S' lying on its side. Of particular interest are the intersection points 
of the two nullclines; these are the fixed points of the dynamics. In the case illustrated in 
Figure 2a, the only intersection is (0,0), an unstable equilibrium. Trajectories intersect 
the s-, resp. a-, nullcline in a direction parallel to the a-, resp. s-, axis. 

The study of the nullclines is of interest because it is often possible to predict how 
a parameter change will affect the dynamics of the system by reasoning about how the 
nullcline diagram will change; the bifurcation we shall be mostly interested in is associ- 
ated with a conspicuous change in this diagram. Note that the s-nullcline is affected by 
parameters w EE and w E1 , whereas the cr-nullcline is affected by parameters w 11 and w 1E . 

(Insert Figure 2 around here) 

Let us consider first the changes brought about by letting parameter w EE grow, starting 
from the point w EE = 12 for which the system oscillates; other parameters are unchanged. 
When w EE grows, the slope of the central, quasi-linear, part of the s-nullcline increases 
(see equation |); that part of the curve rotates about the symmetry center (0,0). As a 
result, the peak of the s-nullcline to the right approaches the upper part of the sigmoid- 
shaped cr-nullcline, and the minimum of the s-nullcline to the left approaches the lower 
part of the cr-nullcline. Eventually, at a certain critical value wf E (subscript 'sn' stands 
for 'saddlenode' — see below), the two curves become tangent to each other. This happens 
in two points at once, near the upper right-hand corner and near the lower left-hand 
corner, due to the symmetry of the system. This situation is depicted in Figure 2b: w EE 
is exactly equal to the critical value wf E (with parameters as above, wf E ~ 14.22), and 
the nullclines are just tangent to each other. 

When w EE grows a little further, each point of contact splits into two intersection 
points, of which one is an attractor. Figure 2c shows this situation, with w EE = 15, 
somewhat above the critical value wf E . Four trajectories are shown, in addition to the 
two nullclines. The system has five fixed points, three unstable ones and two stable ones 
(attractors). Only the stable fixed points are of interest to us; they are very near the upper 
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right-hand and lower left-hand corners of the square, corresponding to high, respectively 
low, excitatory and inhibitory activities. 

The bifurcation occurring at wf E is of the saddlenode type. It results in a drastic 
change of behavior of the system: the periodic attractor disappears and is 'siphoned' into 
the two new point attractors. These two points attract the entire square, except a set of 
measure which includes the three unstable fixed points. Thus, altough this bifurcation 
is caused by a mere local change, namely the intersection of the nullclines, it results in a 
reorganization of the dynamics that is both abrupt and global.^ 

Having described the breakdown of oscillations when parameter w EE is increased, 
we now consider the opposite change, that is, we let w EE decrease. This results in a 
decrease of the slope of the central, increasing, portion of the s-nullcline (equation 
Eventually, the curve becomes monotonically decreasing; this does not alter the number 
of intersections of the nullclines, point (0, 0) remaining the sole equilibrium. However, 
the amplitude of the limit cycle decreases along with w EE . The cycle eventually collapses 
to a point; the equilibrium (0, 0) has then become stable. This can be seen in a linear 
stability analysis of system ^ around point (0,0). It is easily shown that, in case there 
are two complex conjugate eigenvalues^ the real part of these eigenvalues is negative if 
and only if w EE < w 11 + 4T. Thus, w 11 + 4T is a critical value for parameter w EE . We 
define w EE { = f w 11 + AT (with the current parameter setting, w EE { = 6). The change of 
behavior occurring at w EE p{ is a norma^ Hopf bifurcation. 

So far, we studied the behavior of system [| for different values of parameter w EE , 
all other parameters being fixed. In other words, we described the system's behavior 
on a particular 1-dimensional subspace of the 4-dimensional parameter space. We now 
extend this study to a 2-dimensional subspace, the (w EE ,w lE ) plane. Figure 3a is the 
bifurcation diagram of system ^| in that plane, with other parameters as before (w El = 10, 
2 As mentioned, two distinct saddlenode bifurcations take place simultaneously. Such a double bifur- 
cation is not generic; it occurs here due to the symmetry that we introduced when reducing system [l] 
into system ||. 

3 The condition for this is Aw E1 w lE > (w EE + w 11 ) 2 . 

4 That is, supercritical. However, for very large values of u> IE , the bifurcation is subcritical — see 
footnote 5. 
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w 11 = 2). This diagram shows three distinct regions, corresponding to three different 
attractor configurations; unstable fixed points and unstable limit cycles are ignored in 
this diagram. In the middle region — which we call region V, for Periodic — the system 
oscillates. The boundary of this region to the right is the saddlenode bifurcation curve, 
which we denote S; as discussed above, the rightmost region has two point attractors, and 
we call it region T. The leftmost region, which we call 0, has only one point attractor, the 
center of symmetry (0, 0); it is separated from region V by the Hopf bifurcation curve, a 
vertical line of equation w EE = w EE p{ . The curve in the lower left of the diagram, separating 
region O from region T, is the locus of a pitchfork bifurcation. This bifurcation diagram, 
obtained for one particular set of values of the parameters w El , w 11 and (3, is representative 
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of the general case.Q 



(Insert Figure 3 around here) 

In sum, the (w EE , w lE ) bifurcation diagram for system |] is characterized by a central 
periodic-attractor region, a large vertical patch extending to +00 in the w lE direction 
(phase V), flanked by point-attractor regions on each side (phases O and T). The tran- 
sition from V to T is abrupt (<S line), while the transition from O to V is smooth. As 
mentioned in the Introduction, system [l] — the full system — is not amenable to such a 

thorough analysis; however, we shall see in Section [5] that the two systems behave in 
J It is however simplified in two ways. First, the transition from region V to region T is of the 
saddlenode type only for large enough values of w IE ; this range of values corresponds roughly to the 
straight portion of curve <S (Figure 3a). To see why this is so, consider again Figure 2b, the nullclinc 
diagram at the bifurcation, with w = 8. Note that the points of contact between the nullclines appear 
near the corners of the square, far from the origin; this is due to the fact that ui IE is relatively large, hence 
the slope of the er-nullcline at the origin is larger than the slope of the s-nullcline. The bifurcation is then 
of the saddlenode type, as described. If however w lE is small, hence so is the slope of the cr-nullcline at 
the origin, the transition from V to T as w EE is increased takes place differently. A pair of intersection 
points between the nullclines first split off from the origin; these are unstable equilibria. As w EE increases, 
these two equilibria move away from the origin, while remaining inside the large stable limit cycle. At a 
certain critical value for w EE they become stable — a (double) subcritical Hopf bifurcation — and almost 
immediately thereafter the large limit cycle disappears. Thus, the transition from region V to region 
T really takes place in two steps, giving rise to a three- attractor behavior: the system has one large 
limit-cycle attractor as well as two point attractors, the latter being inside the cycle. The region of the 
(u; EE ,ui IE ) plane where this behavior takes place is a strip extending along the lower, curved, part of 
the V/T boundary; it is too narrow to be seen in Figure 3a. (With parameters w E ^ and w 11 as above 
and w lE = 2.75, the three-attractor behavior occurs for w EE between 8.993 and 9.030. For some other 
values of w EI and this behavior does not occur at all, and the transition from V to T is always of the 
saddlenode type.) For the purpose of this paper (see footnotes 8 and 11) it is important to note that the 
point attractors appear either exactly or almost at the same time as the periodic attractor disappears. 
The second approximation in the bifurcation diagram, mentioned only for the sake of completeness, 
concerns the O-to-'P transition. This is generally a smooth, supercritical, Hopf bifurcation. However, as 
mentioned in footnote 4, this Hopf bifurcation becomes subcritical for very large values of iu IE . There 
is thus a narrow region to the left of the bifurcation line w EE = w EE ^ where the limit-cycle attractor 
coexists with the point attractor (0,0); for instance, at w lE = 100, the width of this region is w 0.63. 
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much the same way under the plasticity rules that we shall now introduce. 



3 The regulation equations 

Whereas in the previous section the synaptic weights w EE and w 1E were fixed parameters, 
they will now be made to evolve. Their evolution will obey a Hebbian covariance rule, 
hence be a function of second-order temporal averages of the dynamic variables s and 
a. Synaptic plasticity creates a regulation loop: changing the parameters affects the 
dynamics of the system, which in turn alters the second-order moments of s and a. 
Formally, the regulation is implemented by introducing additional differential equations, 
coupled to system |3] (or to system [TJ — see Section |^). The rate of change of w EE and w lE 
will typically be several orders of magnitude slower than that of s and a. 
Let us first define, for any function of time r(t), a moving time average: 

f(t) = p f r(u)e p{u - t] du. (6) 



Parameter p is a positive constant, physically an inverse time; the larger p, the narrower 
the averaging kernel. Equivalently, f(t) may be defined by a differential equation, more 
convenient for simulation purposes: 

df(t) 



dt 



p(r(t)-f(t)). (7) 



Consider now, with reference to the original stochastic model (Section ^), the instanta- 
neous covariance between two excitatory neurons i and j, defined as: c EE (t) = f (xf(t) — 
xf(t))(x E (t) — xf(t)). If we take the population average (c EE (t)) of this instantaneous 
covariance, we obtain, in the thermodynamic limit iV — > oo, the instantaneous variance 
of s(t): 

c EE (t) ^ (s(t) - s(t)f. (8) 

It is this quantity c EE that we use to regulate the excitatory-to-excitatory synaptic 
weight w EE . The regulation equation is linear in c EE : 

=£ EE (c EE (t) _£ EE) _ (g) 
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Parameters e EE and 9 EE are positive. Note that the quantity c EE (t) is always non- negative; 
the term — 9 EE is therefore necessary to allow for decreases of w EE . 

We shall also consider a regulation for w lE , the synaptic weight from excitatory to 
inhibitory neurons, although this regulation will play a less important role than that of 
w EE . The modification rule for w lE has the same form as equation ^, yet it uses the 
excitatory-to- inhibitory instantaneous covariance, defined as: 

c m (t) d ^(s(t)-s(t))Ht)-m- (io) 

The regulation equation for w 1E then reads: 

^& = £ IE ( C IE (t)-0 IE )- (11) 

In equation [II], 9 lE is a positive constant, as 9 EE in equation^. However, the modification 
rate constant e lE is negative. The main reason for this will be given in the next section; 
for now, note that this choice is consistent with the spirit of Hebb's principle, for, when 
considered postsynaptically to the target neuron, the effect of synapse reinforcement if 
that target neuron is inhibitory is the opposite of the effect obtained if the target neuron 
is excitatory. 



4 Behavior of the regulated reduced system 

This section describes the behavior of the regulated reduced system. We demonstrate that 
each of the two regulation loops introduced in Section |3], when acting separately, brings 
the system to the critical surface S, the locus of an abrupt phase transition (saddlenode 
bifurcation). We then examine the behavior of the system with the two regulation loops 
active simultaneously; we show that under some conditions the state converges to a point 
on S with a remarkable nullcline configuration. 

Before we consider the regulation proper, let us examine how the covariances change 
across the (w EE , w 1E ) plane. Figure 3b shows the values of c EE , the time average of the 
instantaneous variance of s(£),[] along several horizontal lines in the (w EE ,w 1E ) plane. As 
6 This corresponds, in the original system, to the population- and time-average of the covariance, 
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expected, c EE is positive only in region V, where the dynamics is periodic;[] although not 
shown, the same is true of c IE , the time average of the E-to-I covariance. Note that as u> EE 
crosses the O-to-V boundary (Hopf bifurcation) from left to right, c EE increases smoothly 
from to positive values: as discussed above, the amplitude of the limit cycle at this 
bifurcation is infinitesimal. In contrast, the change in c EE and in c IE at S (saddlenode 
bifurcation) is a sharp one, as the system undergoes there a transition from a large limit- 
cycle regime to a fixed-point attractor. 

We now start our study of covariance plasticity by regulating parameter w EE in system 
|3] while all other parameters, including w 1E , remain fixed. The system under study then 
consists of coupled equations |3|, §, ||. Equation |9| prescribes an increase of w EE when 
c EE > 6 EE , and a decrease when c EE < 9 EE . Referring to Figure 3b, we see that to 
the left of S, where c EE is high, the first of the two conditions applies; in this region 
w EE increases. To the right of S the covariance vanishes, and w EE decreases. Therefore, 
w EE (t) is attracted to the transition line S.f\ 

(Insert Figure 4 around here) 

The behavior of this w EE regulation loop is illustrated in Figure 4a for the following 
setting of parameters: w E1 = 10, w 11 = 6, p = .1, 9 EE = .01, e EE = .01. This figure focuses 
on a small region of the (w EE ,w lE ) plane, and shows the projection of the trajectory of 
(c EE ); the latter becomes c EE in the thermodynamic limit N — > oo. In the regulation equation, we use 
the instantaneous covariance c EE (t) rather than its time average c EE (see Discussion). The time-averaged 
variance c EE is used here for illustration purposes only. In order to obtain an essentially constant value 
for c EE rather than an oscillating function of time, different values of p are used for the two averaging 
operations: the kernel used to average c EE into c EE is ten times broader than the kernel used to compute 
s from s. 

7 In general, positive average covariance across a neuronal population indicates collective fluctuations; 
in our simplified two-dimensional system, the only possible nontrivial asymptotic behavior is periodic 
oscillation. 

8 The control parameter 9 EE should be smaller than the value of c EE immediately to the left of S. 
The portion of the boundary line where the bifurcation is a subcritical Hopf rather than a saddlenode 
(footnote 5) yields similar behavior, since the disruption of the large-amplitude limit cycle occurs very 
near the emergence of point attractors (see also footnote 11). 
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(s, a, w EE , w 1E ). Several trajectories are shown, all horizontal since w 1E is a constant, 
These trajectories terminate on the critical line S, and the behavior of the s and a 
components on them is as follows. On the trajectories coming from the left, in the 
V region, (s, a) moves along a cyclic orbit, whose amplitude grows as w EE increases and 
approaches the bifurcation line. On the trajectories coming from the right, in the T region, 
(s, a) stays in one of the two point attractors while w EE decreases until it reaches the 
bifurcation curve. When S is reached, either from the left or from the right, motion does 
not really stop. Rather, w EE sets in a periodic oscillation of small amplitude synchronized 
with a large-amplitude periodic motion of (s, a); the frequency of this oscillation is several 
orders of magnitude slower than in V, hence covariance is small — it matches, on average, 
the control parameter 9 EE . When in this regime, the system spends a long time in one 
of the two almost-attracting corners of the [—.5, +.5] 2 box before leaving it and moving 
rapidly to the other corner. This results in an almost-square wave, a behavior that is 
intermediate between the fast periodic motion observed in V and the bistable situation 
prevailing in T. The period of this oscillation and the amplitude of the oscillation of w EE 
depend on parameters p, 8 EE , and e EE .f\ 

We next consider the w IE -regulated system, where w EE and all other parameters remain 



fixed. This system consists of coupled equations y, [10], |TT|. As noted, the E-to-I covariance 
c IE vanishes outside region V, just like c EE ; within V it varies, in a first approximation, 
like c EE . Since we chose e 1E to be negative, w IE decreases in V and increases in T, whereas 
the opposite was true of w EE when it was regulated. Figure 4b shows this w lE dynamics 
in the same region of the (w EE ,w 1E ) plane as before. Parameters are w El = 10, w 11 = 6, 
p = .1, 9 m = .01 and e lE = —.01. The trajectories are now parallel to the w lE axis, and 
(w EE , w IE ) is again attracted to the critical line S separating region V from region T. 
This is true only to the left of the vertical asymptote of that curve; trajectories to the 
right of that line go to +oo. 

In sum, regulation of either one of the two parameters w EE , w 1E has the effect of 
9 Not shown on Figure 4a is the leftmost part of region V, near the Hopf bifurcation, where the 
limit cycle is of small amplitude hence the condition c EE > 9 EE is not realized. When initialized there, 
the system does not converge to S. However, in both the w EE and the w lE directions, the domain of 
attraction of S extends to +oo. 
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bringing the system to the critical surface S separating the region of oscillation from the 
region of bistable steady firing; the nullcline diagram is then as in Figure 2b. Note that 
when the system is on S, a small perturbation in the weights will elicit either oscillation, 
constant firing at near-maximum rate, or constant firing at near-minimum rate. 

We now turn to the behavior of the system when the two regulation loops act si- 
multaneously; we thus study the system of coupled equations [| ||, [TU], [IT]. Figure 4c 
shows the (w EE ,w lE ) dynamics for the same parameters as before, i.e., w E1 = 10, w 11 = 6, 
p = .1, 9 EE = .01, e EE = .01, 8 m = .01 and e 1E = —.01. It appears from this diagram 
that the evolution proceeds in two clearly distinct stages. In the first stage, which could 
be predicted from the study of the regulation loops acting separately, (w EE , w IE ) moves 
toward line 5.0 When this line is reached, motion slows down considerably — typically by 
several orders of magnitude — and proceeds along the critical line, eventually converging 
to a point on S denoted G in Figure 4c. As before, attractor G is in reality a slow limit 
cycle, of small amplitude in w EE and w lE , and large amplitude in s and a. All four vari- 
ables, s(t), <j(t), w EE (t), and w lE (t), are now synchronized; the distinction between slow 
and fast variables has thus vanished. The basin of attraction of G in the (w EE , w 1E ) plane 
roughly consists of the union of the two domains of attraction of S for the separate w EE 
and w lE regulation dynamics; only the region to the left of and around the Hopf line is 
not attracted to the saddlenode line S and eventually to G. 

(Insert Figure 5 around here) 

The location of S in the (w EE , w IE ) plane depends on the values of the fixed parameters 
w El and w 11 . The location of the attractor G on S further depends on the control param- 
eters 9 EE and 9 lE . When the latter are given identical values, as in the case illustrated 
in Figure 4c, the attractor G has the remarkable property that the s- and a-nullclines 
stand in near overlap over a large portion of the interval [-.5, +.5] (Figure 5); the flow of 
the system in this configuration nearly vanishes on a large one-dimensional manifold in 
10 The direction of this linear motion is roughly parallel to the line w EE = —w lE . This is because 
e EE = — £ IE , 9 EE = 9 lE , and the two covariances c EE and c IE are nearly the same. Another choice of 
parameters would result in a different slope, but otherwise similar behavior. 



16 



the two-dimensional phase space. Further, s(t) and a(t) remain nearly identical at all 
times. p 1 ] 

5 Behavior of the regulated full system 

Recall that system 0, which we used so far, was derived from system |I| by eliminating the 
firing thresholds h E and h 1 (equations 0) in such a way as to make (.5, .5) — (0, 0) in system 
|H — a center of symmetry of the dynamics. While easier to analyze, the reduced system 
is less realistic. There is no clear biological justification for this hard-wired symmetry; 
moreover, when the system is in phase T, i.e., to the right of the critical surface S, it can 
stay for arbitrarily long periods of time in one of the two fixed point attractors, e.g. in 
the high-activity one; this is unrealistic. 

In this section we consider a biologically more plausible way of introducing symmetry 
in the dynamics. Rather than eliminating the thresholds according to equations we 
regulate them, thereby implementing a form of 'soft' symmetry. Regulating the firing 
thresholds in a neural network is a simple way to maintain the mean activity around 
an intermediate, useful, value. This may be viewed as a simplification of the regulation 
mechanisms at work in real brains, which, in all likelihood, involve systems of inhibitory 
neurons acting on various time scales. 

The simultaneous regulation of four parameters results in a complex dynamics, which 
makes a thorough analysis impractical. We shall proceed as follows. We first consider, in 
system [I], the regulation of w EE and h E for a given setting of all other parameters. We 
11 Giving different values to parameters 9 EE and 9 lE mostly affects the behavior of the system after 
it has reached S; if 6 EE is larger, resp. smaller, than 9 1E , the state moves downward, resp. upward, 
on S. When {w EE ,w IE ) is on <S but above point G, the nullclines are tangent to each other but do not 
overlap; such a situation is illustrated in Figure 2b. When (w EE ,w lE ) is on S but below point G, the 
nullclines do overlap, but over a smaller domain. With 9 EE = .0118 and 9 IE = .0100, the state stabilizes 
in the narrow three- attractor region described in footnote 5. The state (s, a) then visits each of the three 
'attractors' in turn: its motion consists of a succession of large-amplitude oscillations (periodic attractor) 
and of spiraling orbits around two symmetric points in the interior of the large cycle (point attractors). 
The amplitude of the motion of (w EE ,w^ E ) remains small. This is a mildly chaotic behavior; a more 
pronounced chaotic behavior will be described in the next section for the full system. 
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show that the system converges to the intersection of two critical curves, each of which 
corresponds to the establishement of one point of contact between the nullclines. We next 
consider the system with all four parameters h E , h , w EE and w 1E regulated, and study 
the projection of the dynamics on the (w EE , w 1E ) plane. There are again two stages; the 
first essentially reproduces the behavior observed with the sole (w EE , h E ) regulation, while 
the second is analogous to that observed when regulating w EE and w lE in the reduced 
system; this applies for a broad range of the remaining fixed parameters w EI and it; 11 . 

Figure 6a is the bifurcation diagram of system [L] in the (it; EE , h E ) plane, for the following 
values of the fixed parameters: it; EI = 10, w lE = 10, w 11 — 1, h 1 — 5. As before, we ignore 
unstable equilibria and unstable limit cycles. As before there are three regions, denoted 
respectively by O, T and V , corresponding to three types of asymptotic behavior: single 
fixed-point attractor; two fixed-point attractors (high and low activity); one periodic 
attractor. We now however subdivide region O — somewhat arbitrarily — according to the 
location of the fixed-point attractor in the phase space: the three subregions denoted Oh, 
O m , and Oi, correspond, respectively, to high, middle, and low activity for this attractor. 
The transition between region V and region O m takes place through the familiar, smooth, 
Hopf bifurcation. The transition between V and Oh, as well as its continuation between 
Oi and T, takes place through a saddlenode bifurcation. We denote by Sh the locus 
of this transition; it marks the appearance of a point of contact between the nullclines 
near the high-activity corner, and is thus similar to the S transition in the reduced 
system. However, due to the symmetry of that system, another point of contact appeared 
simultaneously near the low-activity corner, giving rise to a double bifurcation. In system 
[l] this is no more the case, and the intersection of the nullclines near the low-activity 
corner gives rise to a distinct saddlenode bifurcation line, the transition between V and 
0\, which we denote S\. 

(Insert Figure 6 around here) 

When regulating it; EE according to equation |^ and leaving all other parameters fixed, 
the behavior of system [I] is as follows. When starting in region V to the left of the critical 
line Sh, the system oscillates, covariance is high, hence it; EE increases until it reaches the 
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critical line Sh- A point of contact is then established near the high- activity corner of 
the square. The system settles in a slow periodic attractor, of small amlitude in w EE and 
large amplitude in (s, a), whereby nearly all the time is spent in the high-activity state. 

We now regulate the threshold h E as well, in such a way as to stabilize s, the time 
average of s, around a given target value 9 E : 

^ = e E (s(t)-0 E ). (12) 

The rate constant e E is positive and small, and the control parameter 9 E is chosen well 
in the interior of the interval (0, 1), e.g. between .2 and .8 (remember that in system p] 



the activity variables s and a lie in the interval (0, 1)). To see how equation [12| achieves 
the desired regulation, note for instance that, if s(t) > 8 E , h E will increase, which in turn 
will result in a decrease of s(t). 

When both w EE and h E are regulated, the system converges to the intersection of the 
two critical lines Sh and Si. In effect, we saw that the full system, when at a generic point 
of Sh, stays nearly all the time in the high-activity state; this results in a high value of s. 
To achieve the condition s ~ 8 E , the equilibrium for equation VI, the system can only be 
on Si as well. 

The joint (w EE , h E ) dynamics is illustrated in Figure 6b, for parameters w El , w 1E , w 11 
and h l as above, and p = .2, e E = .001, 9 E = .5, e EE = .01, 9 EE = .01. The intersection of 
Sh and Si, denoted F in Figure 6b, is reached from all directions in the (w EE , h E ) plane. 
When coming from low w EE values, the system oscillates and converges to F through 
region V. When coming from high w EE values, the system reaches F through region T, 
where it bounces back and forth between the high- and low-activity point attractors (an 
oscillation much slower than in V). 

The nullcline diagram for point F of Figure 6b is illustrated in Figure 6c. There are 
now two points of contact between the nullclines, a situation more degenerate than the 
one that obtains from regulating w EE only, but 'equivalent' to the situation obtained in 
the reduced system by regulating a single parameter, w EE or w lE (compare Figure 6b to 
Figure 2b). What characterizes the dynamics at point F is that the system is on the 
verge of oscillation and on the boundary of each of the two steady-firing phases. 

We finally consider the system with the four parameters h E , h l , w EE and u> IE regulated. 
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We thus include, in addition to equations [T], H, H [K| |TT] and [12], a regulation equation for 
the inhibitory threshold h 1 : 

^-=s\a(t)-O l ). (13) 

As in equation [12], the rate constant e l is positive and small, and 9 l is chosen in the 
interval (.2, .8), with 9 l « 6* E . The variables now include the activity state (s, a) as well 
as the four regulated parameters h E , h , w EE and w lE . 

Figure 7 illustrates the behavior of this system projected on the (w EE , w 1E ) plane, 
for the following parameter values: w El = 10, w 11 = 6, p = .05, 9 EE = .01, e EE = .01, 
6 lE = .01, e lE = -.005, 6 E = .5, e E = .005, 9 1 = .5, e 1 = .002. In the sequel, this parameter 
setting will be referred to as standard. In a first stage, the system converges to a doubly 
critical point F as described above; each such point F belongs to the common boundary 
of the regions of oscillation, high steady firing, and low steady firing. Although we cannot 
thoroughly characterize the surface of F points in the four-dimensional (w EE , w lE , h E , h l ) 
space as we did in the (w EE ,w IE ) plane for the reduced system, there is, as remarked 
above, a functional equivalence with the S surface. Note that the projection of the F 
surface on the (w EE ,w 1E ) plane has a shape quite similar to that of S in the reduced 
system. As before, when the system reaches a point F, all variables settle in a slow, 
synchronous, almost-periodic motion. The oscillation of s and a is a nearly rectangular 
wave, the system spending nearly all its time in the two corners of the square, where 
the relative amount of time spent in each corner is determined according to the value of 
parameter 9 E (~ l ). As before too, the first stage, which consists of the convergence to 
a doubly critical point F, is robust against parameter changes; most parameters can be 
individually varied over several orders of magnitude without qualitatively affecting this 
part of the behavior. 

The second stage, consisting of a much slower motion on the F surface, depends on 
the values of the various parameters. For most parameter settings, including the standard 
set (see above), the behavior on this critical surface is a slow, simple, periodic motion, of 
large amplitude in (s,a) and very small amplitude in (w EE ,w lE ). The system eventually 
settles in a periodic attractor of this simple type, denoted again G in Figure 7. Figure 8a 
shows the (s, a) projection of this attractor for the standard parameter set; its (w EE , w IE ) 
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projection is a small cycle around point G, whose nullcline diagram is similar to the one 
shown in Figure 5 (largely overlapping nullclines). 

There exists however a small region of parameter space, mostly around e l ~ e E , for 
which a variety of more complex behaviors are observed during the second stage. The 
following two cases are examples of such complex behavior. For parameters as above 
(standard) except that e E = .0051, e l = .0046, and 9 EE = .011, the system settles in 
a complex quasi-periodic motion (Figure 8b). For parameters as standard except that 
e l = e E = .005, the system displays strongly chaotic behavior (Figure 8, c-e).0 Both 
of these behaviors are actually attractors, reached after considerable time, yet similar 
behaviors also take place while the system is still moving slowly on the critical surface. 

To summarize, both in the reduced and in the full system, convergence to a doubly 
critical surface between the regions of fast oscillations and of high and low steady firing 
takes place reliably for a broad range of parameters. Once this doubly critical surface 
is reached, motion becomes slow, depends on parameters, and, when examined in detail, 
reveals a variety of behaviors, ranging from simple periodic firing to chaos. 

6 Discussion 

This paper proposes that a regulation mechanism underlies criticality in brain dynamics. 
In such a scheme, regulation stabilizes the dynamics near an instability. The force driving 
the system towards criticality is a covariance-governed modification of synaptic efficacies in 
a recurrent network. Although it has been argued that criticality in some physical systems 
may be self-organized (Bak et al. 1987), this phenomenon may not be very widespread. 
The nervous system is moreover regulated homeostatically to withstand perturbations 
of various sorts. It is then of interest to explain how criticality in brain dynamics may 
nevertheless arise, from a specific, well-documented, mechanism of synaptic plasticity. 

The chief motivation for viewing the dynamics of the nervous sytem as critical is the 
observation that brains are very sensitive organs. Not only do our brains draw distinctions 
12 This behavior takes place only for some initial values in the (w EE , w IE ) plane; other initial values 
converge to a point attractor. 
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between stimuli that differ only in minute details, they also allow us to establish subtle 
yet clear-cut boundaries between cognitive categories at different levels. The various 
manifestations of hyperacuity in sensory systems may be no more than elaborate forms 
of signal amplification; however, in higher cognitive functions such as language, abstract 
thinking, or, for instance, artistic composition, the ability to create a new category by 
drawing a fine line — an ability which manifests itself early in life and stays with us for 
a long time — argues in favor of regulated criticality. Such a mechanism appears to be 
necessary in order to explain how sensitivity is maintained in the face of the profound 
changes that affect the connectivity of the brain throughout development and learning. If 
no such mechanism were present, one would expect that the ongoing modification of the 
networks which carry mental representations would soon bring these networks to generic 
states; as mentioned in Section |IJ, a dynamical system in a generic state does not show 
high susceptibility to external influences. 

The emergence of new cognitive categories may in effect be likened to a process of 
morphogenesis in embryology, or differentiation in cell biology. A biological structure that 
is about to undergo differentiation is at that particular instant of time unstable, and, as the 
well-studied mechanism of induction shows, highly susceptible to external signals. From 
a dynamical-system perspective, the emergence of qualitatively new behavior, e.g., the 
splitting of one attractor into two, is a bifurcation. The complexification of an individual's 
cognitive apparatus in the course of his or her life may be viewed as an open-ended 
sequence of such bifurcations. Such an interpretation has been defended by Rene Thorn 
(1975), and related ideas have been expressed by several authors (e.g. van der Maas 
and Molenaar 1992). Thorn (1975) also suggested that structurally stable non-generic 
singularities may arise from a process he termed the stabilization of thresholds; this process 
itself would result, in various biological contexts, from the reinforcment of homeostatic 
mechanisms.0 

The covariance plasticity rule we use is linear and straightforward. Equation |5] may 

be viewed as a mean-field version of the covariance rule used in the associative-memory 
13 We thank Jean Petitot for pointing out to us that regulated criticality as proposed here is closely 
related to Thorn's ideas. 
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literature (see e.g. Willshaw and Dayan 1990). However, we make a rather different use of 
this rule. In an associative-memory model, pre- and post-synaptic activities are generally 
assumed to be independent, yielding a zero expected value for the covariance. Weights are 
modified according to the instantaneous covariance, and, as noted in Dayan and Sejnowski 
(1993), storage is marked by the departures of the empirical average of this quantity from 
its expected value, which is zero. In our model, the expected covariance is positive in the 
oscillatory phase. The regulation mechanism acts on a slow time scale, and, although we 
use the instantaneous covariance in the modification rule (Equation [|), we might as well 
have used the time-averaged covariance; fast variations of the instantaneous covariance 
are actually smoothed out in the integration of the differential equation. Of course by the 
very principle of regulation proposed, the system does not dwell in the oscillatory phase; 
in the regulated state, the average covariance is low. 

The other major difference between the situation studied here and the associative- 
memory paradigm is the assumption of uniform weights. As noted in Section ^, the co- 
variance in our uniform-weight network is simply the variance of the population-averaged 
activity about its mean, and is always non-negative. This makes it necessary to subtract 
from it a positive constant 8 EE in order to allow for decreases of the weights. Thus, 
whereas in associative- memory models a synaptic weight decreases as a result of negative 
instantaneous covariance between the pre- and post-synaptic neurons, the condition for 
weight decrease in our model is that the mean covariance be small or zero, which happens 
when the system is at rest in a point attractor, of either low or high activity. 

The uniform-weight network used in the present study lends itself to a detailed math- 
ematical/numerical analysis. We performed a bifurcation analysis of the continuous-time 
differential system that describes the behavior of this network in the thermodynamic limit. 
This analysis (Section [|) reveals, among other features, the existence of a critical surface 
Sh in parameter space, where the system undergoes an abrupt transition from oscillatory 
behavior to high-rate steady firing. We showed (Section ^) that Hebbian modification of 
the E-to-E synaptic weights drives the system toward this surface Sh', this is the main 
mechanism of regulated criticality proposed. 

However, when the system is at a general position on Sh, it spends most of its time 
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in the high-activity state, which is undesirable. A more realistic situation results from 
introducing some form of symmetry between the high- and low-activity phases, making 
use of the firing thresholds. We investigated two ways — formally different but function- 
ally equivalent — to do this. The mathematically simpler way is to enforce an accurate 
symmetry on the dynamics, by imposing an appropriate relationship between the firing 
thresholds and the synaptic weights (equations 0). This results in a reduced system, with 
only four parameters; in this system, there occurs a double bifurcation when the system 
traverses the critical surface — now denoted S — separating the oscillatory phase from the 
bistable, high/low, steady-firing phase (Section g). Regulation of the sole E-to-E weight 
brings the system to this doubly critical surface S (Section ^). 

A biologically more satisfactory solution is to regulate one or both of the firing thresh- 
olds so as to control the mean firing rates (Section Thus, when we regulate the thresh- 
old for the excitatory neurons in addition to the E-to-E weight, the system converges to 
the intersection of Sh with another critical surface, Si, which separates the oscillatory 
phase from the low-activity fixed-point region. Intersection points between Sh and Si are 
again doubly critical, and they attract the system for a wide range of parameter values. 

When the system is on this doubly critical surface, it takes only a small weight per- 
turbation to induce either of the three behaviors: intrinsic oscillation (region V of section 
d), high activity, quiescence. It is easily seen that, when in this state, the network can 
also be efficiently driven by a small-amplitude time- varying signal, i.e., an external field; 
it is thus highly sensitive to input. The stabilization at the boundary of a region of oscil- 
latory behavior appears to be consistent with at least one conclusion that can be safely 
drawn from the recent literature on cortical oscillations, namely the fact that the precise 
conditions under which these oscillations occur are difficult to pin down. 

We further investigated the effect of regulating the E-to-I weight in addition to the 
E-to-E weight, according to a similar covariance rule. We showed that regulating these 
two weights as well as the two firing thresholds results, under appropriate parametric 
conditions, in convergence to an even more degenerate state. When the system is in that 
state, its flow vanishes on an entire one-dimensional curve in the two-dimensional phase 
space, instead of on isolated points. This convergence is slow and parameter-dependent, 
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yet it is interesting to note that when the system is in or near this highly degenerate state 
it exhibits a range of diverse behaviors, including chaos (Section [5]). The chaotic behavior 
shown in Figure 8, c-e consists of an irregular sequence of spontaneous transitions between 
the three fundamental phases of the system: oscillatory, high-activity, low-activity. 

While the uniform-weight network studied in this paper lends itself to a convenient 
mathematical analysis, it would be interesting to know whether critical behavior may arise 
from local covariance plasticity, where synaptic changes are made to depend on pre- and 
post-synaptic activities relative to individual synapses. This question should be focused 
by considerations about the elaborate forms of input sensitivity that could play a role in 
higher brain functions. 
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Figure 1: Mean activities of excitatory and inhibitory populations in a moderate-size 
uniform- weight system exhibiting oscillatory behavior (N = 70; Glauber dynamics). 
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Figure 2: Behavior of mean-field system for different values of the E-to-E synaptic 
weight w EE . Diagrams show trajectories and nullclines. (a) w EE = 12 (parameters 
are as in Figure 1); all trajectories converge to a limit cycle, (b) w EE = w EE ~ 14.22; 
the system is at the saddlenode bifurcation: nullclines are tangent to each other (no 
trajectories shown), (c) w EE = 15; nullclines intersect, the periodic attractor has 
vanished, two point attractors have appeared. 
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Figure 3: Bifurcation diagram in (w EE , w lE ) plane, (a) Diagram shows three regions, 
characterizing different attractor configurations. Region O: single point attractor, 
of intermediate activity level; region V: periodic attractor, as depicted in Figure 
2a; region T: two point attractors, of high and low activity, as depicted in Figure 
2c. Transitions between regions occur through bifurcations, of Hopf type between 
O and V, of saddlenode type between V and T (curve S), and of pitchfork type 
between O and T. (b) Average covariance along four different lines of constant w lE 
in the (w EE ,w 1E ) plane. Note the sharp variation of the covariance on the critical 
line <S separating V from T. 
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Figure 4: Regulation of system |3] by covariance plasticity, (a) w EE is regulated, 
w lE is constant: state converges to critical surface S. (b) w 1E is regulated, w EE is 
constant: state converges to critical surface S. (c) both w EE and w lE are regulated: 
state converges to a particular point, denoted G, on critical surface S. 
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Figure 5: Nullcline diagram at point G (see figure 4c). Nullclines overlap almost 
perfectly over much of the interval [—.5, .5]. 
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Figure 6: Regulation of two parameters in system |l| (a) Bifurcation diagram in 
(u> EE , h E ) plane, (b) Regulation of w EE and h E causes convergence to point F, the 
intersection of critical lines Sh and S[. (c) Nullcline diagram at F. 
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Figure 7: Behavior of full system under simultaneous regulation of four parameters. 
Diagram shows projection on (w EE , w lE ) plane, illustrating the similarity of behavior 
with reduced system (compare with Figure 4c, but note difference of scales). Limits 
of the attraction basin to the left are roughly indicated by the starting points of the 
trajectories shown; attraction basin is unbounded in all other directions. 
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Figure 8: Various behaviors of regulated full system after it has reached critical 
surface (Figure 7). Diagrams show s(t) for three slightly different parameter settings 
(see text); in all cases, the projection of the motion on the (u> EE , w lE ) plane remains 
of small amplitude, (a) Simple periodic attractor, point G of Figure 7; similar 
periodic attractors are reached for most parameter settings, (b) Complex quasi- 
periodic attractor. (c-e) Chaotic attractor; for a given parameter setting, three 
diagrams corresponding to different instants of time and different lengths of time; 
characteristic are the irregular transitions between the high-activity, low-activity, 
and oscillatory phases. 
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